Observation processes

• in practice, often on some sub-areas of the area of interest can be surveyed
• interest in entire observation area

Overview of Distance Sampling

• Distance sampling is a family of related methods for estimating the abundance and spatial distribution of wildlife populations.

• Distance sampling is based on the idea that animals further away from observers are harder to detect than animals that are nearer.

Density surface models

Coupling distance sampling data with spatial modelling allows us to map spatially varying density

• Traditionally, this is achieved in a two-stage approach by (i) using a detectability point estimate to create an offset vector to (ii) use within GLM or GAM models for count response data.

Density surface models

Coupling distance sampling data with spatial modelling allows maps of spatially varying density to be produced.

• Traditionally, this is achieved in a two-stages approach by (i) using a detectability point estimates to create an offset vector to (ii) use within GLM or GAM for count response data.

• This requires binning the data into counts based on some discretisation of space.

Density surface models

Coupling distance sampling data with spatial modelling allows maps of spatially varying density to be produced.

• Traditionally, this is achieved in a two-stages approach by (i) using a detectability point estimates to create an offset vector to (ii) use within GLM or GAM for count response data.

• A major downside to this approach is the propagation of uncertainty from the detection model to the second-stage spatial model.

Thinned Point Process

• flexible approach that can include spatial covariates to model the mean intensity and a spatially structured random effect to account for unexplained heterogeneity not captured by the covariates.

• to account for the imperfect detection of points we specify a thinning probability function \[ g(s) = \text{Prob}(\text{a point at s is detected}|\text{a point is at s}) \]

• key property – an LGCP with intensity \(\lambda(s)\) that thinned by a probability function \(g(s)\), follows also a LGCP with intensity:

\[ \underbrace{\tilde{\lambda}(s)}_{\text{observed process}} = \underbrace{\lambda(s)}_{\text{true process}} \times \underbrace{g(s)}_{\text{thinning probability}} \]

Thinned Point Process

Lets visualize this on 1D: Intensity function with points

Thinned Point Process

Intensity (density) function with points and transect locations

Thinned Point Process

• Detection function \(\color{red}{g(s)}\)

• Here \(\color{red}{g(s) =1}\) on the transects (at x = 10,30 and 50).

Thinned Point Process

• Detection function \(\color{red}{g(s)}\) and detected points

Thinned Point Process

Thinned Point Process

The detection function describes the probability \(\color{red}{p(s)}\) that an point is detected

Thinned Point Process

Thinned Point Process

Observations are from a thinned Poisson process with intensity \(\lambda(s) \color{red}{p(s)}\)

Detection Function

• Standard distance sampling approaches specify \(g(s)\) as a function that declines with increasing distance
  • horizontal distance to the observer for point transects
  • perpendicular distance to the transect line for line transects

Detection Function

• Standard distance sampling approaches specify \(g(s)\) as a function that declines with increasing distance
  • horizontal distance to the observer for point transects
  • perpendicular distance to the transect line for line transects

Detection Function

• Standard distance sampling approaches specify \(g(s)\) as a function that declines with increasing distance
  • horizontal distance to the observer for point transects
  • perpendicular distance to the transect line for line transects
• The thinning probability function is specified as a parametric family of functions.

Detection Function

• Standard distance sampling approaches specify \(g(s)\) as a function that declines with increasing distance
  • horizontal distance to the observer for point transects
  • perpendicular distance to the transect line for line transects
• The thinning probability function is specified as a parametric family of functions.
• The thinned-LGCP likelihood is given by:

\[ \pi(\mathbf{s_1},\ldots,\mathbf{s_m}) = \exp\left( |\Omega| - \int_{\mathbf{s}\in\Omega}\lambda(s)g(s)\text{d}s \right) \prod_{i=1}^m \lambda(\mathbf{s}_i)g(\mathbf{s}_i) \]

• To make \(g(s)\) and \(\lambda(s)\) identifiable, we assume intensity is constant with respect to distance from the observer.

  • In practice this means we assume animals are uniformly distributed with respect to distance from the line

Putting all the pieces together

Putting all the pieces together

Putting all the pieces together

Putting all the pieces together

• The pdf of detected distances is \(\pi(d_1,\ldots,d_m|m) \propto \prod_{i=1}^m\dfrac{\tilde{\lambda}(d_i)}{\int_0^W \tilde{\lambda}(d)\text{d}d}\)

Putting all the pieces together

• The pdf of detected distances is \(\pi(d_1,\ldots,d_m|m) \propto \prod_{i=1}^m\dfrac{ g(d_i)}{\int_o^W g(d) \text{d}d}\) if \(\color{red}{\tilde{\lambda}(d_i)} = \lambda \color{red}{g(d_i)}\)

An approximation: Strips as lines

The intensity at a point \(\mathbf{s}\) on the line becomes \(2W\lambda(s)\) instead of \(\lambda(s)\).

• Let \(\pi(d)\) be the probability that the point is at a distance \(d\) from the line.

• Let \(p(d)\) be the probability that is detected given it is at \(d\).

Then, the thinning probability becomes \(\pi(d)\times p(d)\), assuming the points are uniformly distributed within the strip then \(\pi(d) = 1/W\) (the density of distances is assumed to be constant on the interval \([0,W]\)).

This updates our thinning intensity to

\[ \log \tilde{\lambda}(s) = \underbrace{\mathbf{x}'\beta + \xi(s)}_{\log \lambda(s)} + \log p(d) + \log \times(2/2W) \]

• Typically \(p(d)\) is a non-linear function, that is where inlabru can help via a Fixed point iteration scheme (further details available in this vignette)

Example: Dolphins in the Gulf of Mexico

In the next example, we will explore data from a combination of several NOAA shipboard surveys conducted on pan-tropical spotted dolphins in the Gulf of Mexico.

• A total of 47 observations of groups of dolphins were detected. The group size was recorded, as well as the Beaufort sea state at the time of the observation.

• Transect width is 16 km, i.e. maximal detection distance 8 km (transect half-width 8 km).

Step 1: Define the SPDE representation: The mesh

First, we need to create the mesh used to approximate the random field. We can either:

1. Create a nonconvex extension of the points using the fm_mesh_2d and fm_nonconvex_hull functions from the fmesher package:

Step 1: Define the SPDE representation: The mesh

First, we need to create the mesh used to approximate the random field. We can either:

2. Use a pre-defined sf boundary and specify this directly onto the mesh construction via the fm_mesh_2d function

Step 2: Define the Detection function

We start by plotting the distances and histogram of frequencies in distance intervals.

Then, we need to define a half-normal detection probability function. This must take distance as its first argument and the linear predictor of the sigma parameter as its second:

Example: Dolphins in the Gulf of Mexico

# define model component
cmp = ~ Intercept(1) + 
  space(main = geometry, model = spde_model) +
  sigma(1,
    prec.linear = 1,
    marginal = bru_mapper_marginal(qexp, pexp, dexp, rate = 1 / 8)
  ) 

# define model predictor
eta  = geometry + distance ~ space +
  log(hn(distance, sigma)) +
  Intercept + log(2)

# build the observation model
lik = bru_obs("cp",
              formula = eta,
              data = mexdolphin$points,
              ips = ips)

# fit the model
fit = bru(cmp, lik)

# build integration scheme
distance_domain <-  fm_mesh_1d(seq(0, 8,
                              length.out = 30))
ips = fm_int(list(geometry = mesh,
                  distance = distance_domain),
             samplers = mexdolphin$samplers)

Example: Dolphins in the Gulf of Mexico

# define model component
cmp = ~ Intercept(1) + 
  space(main = geometry, model = spde_model) +
  sigma(1,
    prec.linear = 1,
    marginal = bru_mapper_marginal(qexp, pexp, dexp, rate = 1 / 8)
  ) 

# define model predictor
eta  = geometry + distance ~ space +
  log(hn(distance, sigma)) +
  Intercept + log(2)

# build the observation model
lik = bru_obs("cp",
              formula = eta,
              data = mexdolphin$points,
              ips = ips)

# fit the model
fit = bru(cmp, lik)

# build integration scheme
distance_domain <-  fm_mesh_1d(seq(0, 8,
                              length.out = 30))
ips = fm_int(list(geometry = mesh,
                  distance = distance_domain),
             samplers = mexdolphin$samplers)

Example: Dolphins in the Gulf of Mexico

# define model component
cmp = ~ Intercept(1) + 
  space(main = geometry, model = spde_model) +
  sigma(1,
    prec.linear = 1,
    marginal = bru_mapper_marginal(qexp, pexp, dexp, rate = 1 / 8)
  ) 

# define model predictor
eta  = geometry + distance ~ space +
  log(hn(distance, sigma)) +
  Intercept + log(2)

# build the observation model
lik = bru_obs("cp",
              formula = eta,
              data = mexdolphin$points,
              ips = ips)

# fit the model
fit = bru(cmp, lik)

# build integration scheme
distance_domain <-  fm_mesh_1d(seq(0, 8,
                              length.out = 30))
ips = fm_int(list(geometry = mesh,
                  distance = distance_domain),
             samplers = mexdolphin$samplers)

Example: Dolphins in the Gulf of Mexico

# define model component
cmp = ~ Intercept(1) + 
  space(main = geometry, model = spde_model) +
  sigma(1,
    prec.linear = 1,
    marginal = bru_mapper_marginal(qexp, pexp, dexp, rate = 1 / 8)
  ) 

# define model predictor
eta  = geometry + distance ~ space +
  log(hn(distance, sigma)) +
  Intercept + log(2)

# build the observation model
lik = bru_obs("cp",
              formula = eta,
              data = mexdolphin$points,
              ips = ips)

# fit the model
fit = bru(cmp, lik)

# build integration scheme
distance_domain <-  fm_mesh_1d(seq(0, 8,
                              length.out = 30))
ips = fm_int(list(geometry = mesh,
                  distance = distance_domain),
             samplers = mexdolphin$samplers)

Results: posterior summaries

Results: posterior summaries

Results: predicted densities

Results: predicted densities

distdf <- data.frame(distance = seq(0, 8, length.out = 100))
dfun <- predict(fit, distdf, ~ hn(distance, sigma))

Results: predicted expected counts

We can look at the posterior for the mean expected number of dolphins. Remember, the number of dolphins over the whole domain \(\Omega\) is \[ N(\Omega)\sim\text{Poisson}(E_{\Omega}), \text{ with } E_{\Omega} = \int_{\Omega}\lambda(s)ds \] so the expected mean number of dolphins is \[ E_{\Omega} = \int_{\Omega}\lambda(s)ds = \int_{\Omega}\exp\{\beta_0+\omega(s)\}ds \]

predpts <- fm_int(mexdolphin$mesh, mexdolphin$ppoly)
Lambda <- generate(fit, predpts, ~ sum(weight * exp(space + Intercept)),
                   n.samples = 3000)

Estimate:

Results: predicted expected counts

If want to predict the expected counts We can also get Monte Carlo samples for the expected number of dolphins as follows:

[1] 601  10